从组学数据中进行机器学习

机器学习在组学数据上的应用方向之一是药物发现。

本文为《Artificial Intelligence in Drug Design》一书第18章的内容，本章介绍了如何训练和分析基于组学数据的ML模型，以及一个建立预测药物诱发肝损伤模型的案例。

该书的总体介绍见Springer推出新书《Artificial Intelligence in Drug Design》

摘要

机器学习（ML）已经加速了许多科学领域的发现，是一些新产品背后的驱动力。最近，不断增长的样本量使得ML方法能够在更大的组学研究中得到应用。本章提供了一个如何使用ML对组学数据集进行典型分析的指南。同时，本章展示了一个如何根据转录组学数据（来自LINCS L1000数据集）建立一个预测药物诱发肝损伤模型的案例，涵盖了从数据探索和模型训练（包括超参数搜索）到最终模型的验证和分析的最佳实践和陷阱。重现结果的代码可在https://github.com/Evotec-Bioinformatics/ml-from-omics 获得。

1 简介

在许多机器学习领域，数据是无处不在的。文本和图像可以在网络上找到，顾客在网上购物时也会留下详细的痕迹。相比之下，组学数据是很难得到的，因为它的创建需要高度熟练的工作人员和配备昂贵设备的实验室。

由于本章的重点是机器学习的应用，因此对各种组学技术不作过多介绍。

给定一定数量的生物样本N，用一个组学技术测定一定数量的分子实体P。任何熟悉机器学习的人都会立即认出这是一个特征矩阵M ∈ N × P。通常，P的特征在几百到几千个的范围内。

转录组学领域的一个例子是通过RNA测序（RNA-seq）获得的基因计数矩阵M。每个条目 Mi,j 代表了一个样本 i 的与基因 j 相关的mRNA分子的数量。

有监督的机器学习方法在组学图谱数据上有大量可能的应用。在精准医疗中，从已诊断的病人身上获得的样本可以作为训练样本，以便在未来为病人做出更好的诊断。这对那些难以诊断或治疗昂贵的疾病来说特别有价值。另一个应用领域是药物发现。给出一组经典的化合物及其组学图谱作为训练集，ML算法可以学习预测新的化合物的作用模式，或为新药的设计提供信息。这样就有可能在研究的早期阶段识别潜在的有毒化合物，以避免不必要的体内研究，包括临床试验。

本章使用最近发表的LINCS L1000数据集的一个子集，来说明组学数据集的特点和陷阱。这个数据集的最初目的是提供大量的转录组图谱，阐释不同化合物对各种人类细胞类型的影响。

本章的重点是许多化合物的一个特别令人不快的副作用：肝脏损害。在药物发现中，识别引起药物性肝损伤（Drug-Induced Liver Injury，DILI）的化合物是非常重要的。如果不安全的候选化合物能够在研究过程的早期就被分类出来，就可以减少资源的浪费，降低临床试验的风险。

上述数据集包含了640种不同化合物的6000个转录组图谱，分布于不同的细胞系、剂量和治疗时间。3568个样本被归类为DILI阳性，2432个为DILI阴性。由于对所有样本进行全RNA测序的成本过高，作者选择了978个有代表性的标志性基因进行定量分析。这种有针对性的方法（称为L1000）要便宜得多，因此可以对成千上万的样本进行分析。每个分析只包含那些标志性基因的表达水平。更确切地说，所有的例子都是利用第5级数据（moderated Z-scores）作为机器学习模型的特征。为了具有充分的可比性，将数据集分成4800个训练样本和1200个验证样本。

当前的任务是建立一个二元分类模型，将一个新的化合物转录组特征分类为DILI阳性或阴性。然而，在深入研究模型选择的细节之前，让我们先仔细看看数据。

2 数据探索

作为机器学习的从业者，我们最好对我们的输入数据进行基本的质量检查。至少，我建议计算每个特征的平均数、方差、最小值和最大值等汇总统计数据，并确定缺失值的比例。由于它们不太可能对我们的分类任务有帮助，我们可以立即放弃方差很低或有很多缺失值的特征。此外，我们可以绘制每个特征的直方图，以评估其数值的分布。虽然对于组学数据来说，直方图的数量可能很大，但我们至少可以快速浏览一下，发现任何引人注目的怪异现象，如不同的尺度或严重倾斜的分布。使用这些方法，我们能够熟悉数据并获得一些直觉，从而在模型构建步骤中为我们的决策提供参考。

大量的潜在因素会导致异常值或影响一个或多个样品的系统偏差的引入。这些可能很简单，就像计量错误或实验室的泵出了问题。不幸的是，也有许多微妙的原因。也许在假期中处理最后几个样品的实验室技术人员换成了其他技术人员，遵循了稍微不同的协议，或者实验室的空调在炎热的日子里出现了故障。所有这些因素都会影响测量的结果，从而导致异常值和批次效应，需要去除或纠正。

识别最明显的离群值和批次效应的方法之一是使用降维方法将数据投射到二维，并检查散点图。图1显示了用UMAP在densMAP模式下生成的例子数据的二维嵌入。最引人注目的观察结果是有两个集群：左手边的异质性集群和右手边的以DILI负数为主的集群。简单地说，在我们的分类任务中，有一些容易的案例，也有一些困难的案例。我只能猜测这两个集群是否是用于从原始数据集中选择样本的方法的一个伪装，或者是否有生物学上的解释。至少与任何化合物、剂量、治疗时间或细胞系没有对应关系。一个小得多的集群大致位于坐标（0, 10），主要由用药物Vorinostat治疗的样本组成。总之，数据中没有明显的离群值或批次效应。

图1使用UMAP生成的数据集的二维嵌入

每个点代表一个样本，其颜色表示DILI等级

如果一个数据集包含复制，另一个明显的排除一些错误的方法是检查复制是否相互一致。从视觉上看，这可以在PCA降维图中通过给所有的复制体以相同的颜色来完成。期待的是样本按复制组形成密集的聚类。然而，这种方法并不能适用于大量的样本。在这种情况下，我建议计算复制体之间的相关关系。通常情况下，技术复制体之间的皮尔逊相关度在0.98或以上。一个明显较低的值最有可能表明任何一个样品有问题。

接下来，我们检查特征之间以及每个特征与目标之间的关系。通常情况下，从组学实验中获得的特征会形成相关的集群。同一途径的基因往往是核心关联的，一种代谢物的丰度取决于其前体的可用性。这些强烈相关的特征没有增加什么信息，在估计线性模型的系数时可能会出现问题，并使模型的解释变得复杂。因此，最好是把它们合并成一个特征（例如，通过平均），或者选择一个代表，删除其他的。图2左边的直方图显示了所有特征的最大绝对Spearman相关度的分布。换句话说，我们对每个特征进行检查，哪个特征与其他特征的相关性最大，并在直方图中报告这一绝对相关性。然而，各特征之间只有适度的相关性（<0.8），因此没有必要删除任何特征。这很可能是由于标记基因的选择首先是为了捕捉正交信息。

图2

左图：任何一对特征之间的最大相关度直方图。右图：每个特征与目标类之间的相互信息直方图

使用单变量方法来研究一个特征对于目标类是否提供信息也是值得的。这给出了模型需要多少特征才能实现合理性能的第一个提示。也许有几个高分的特征，可以很容易地进行分类，或者信息被分散到许多特征中。在最坏的情况下，数据中根本就没有可用的信息。然而，具有低单变量分数的两个特征的组合仍然可能对分类任务非常有用。图2的右边部分显示了一个直方图，总结了每个特征和目标之间的相互信息（mutual information, MI），这里的目标是DILI类。

研究选择了相互信息作为单变量测量，因为它也能捕捉到部分和非线性关系。缺点是，它的计算量相当大，因此很可能不适合大数据集。虽然所有的分数都远远低于0.67的最大可能分数（目标向量与自身的MI），但显然有许多无信息的特征和极少数有信息的特征。这意味着某种特征选择应该是机器学习模型的一部分。

3 模型的定义

在数据探索步骤中收集到的知识基础上，我建议使用一个使用Python 3和scikit-learn实现的简短管道作为机器学习模型。它由三个步骤组成。(1）特征标准化，（2）特征选择，（3）支持向量机（SVM）作为分类器。

在当前深度学习大潮的背景下，这种方法可能显得不合时宜。然而，这些老式的机器学习算法的表现与深度神经网络（DNN）的公布结果相当。此外，它们的训练速度更快，更容易解释。特征选择步骤的动机是，在数据探索过程中，我们只观察到少数具有合理单变量得分的变量。SVM假设所有的特征都像标准高斯一样大致分布，均值为零，方差为单位。但对于例子的数据来说，情况并非如此。因此，管道的第一步是对所有特征进行标准化。当然，该管道有一些超参数，包括标准化和特征选择的实际方法、选择的特征数量、SVM惩罚项中使用的规范以及正则化参数的值。

正如介绍中所指出的，获得大量的样本是相当昂贵的，因此，在处理组学数据时，我们经常处于N≪P的状态。这里存在着声名狼藉的"维度诅咒"。在最坏的情况下，每个样本都有一个独特的特征，而模型只是学习使用这些特征来列举它以前见过的样本。很明显，这导致了很差的概括性，因此在独立测试集上的表现将是灾难性的。因此，我们需要格外小心，不要让我们的机器学习模型过度拟合。正则化和特征选择是实现这一目标的最常用工具。因此，在超参数搜索中，要选择的特征数量和正则化强度参数都将被涵盖。

4 超参数搜索

在这里，我使用一个简单的网格搜索来寻找超参数的最佳值。由于模型的训练速度很快，而且超参数的数量也不多，这在计算上还是可行的。然而，具有许多超参数的复杂模型可能需要更复杂的搜索策略，如Hyperband或Bayesian方法。作为特征选择方法，该管道使用相互信息（如数据探索步骤）。此外，我将SVM的最大迭代次数增加到50,000次，以确保收敛性。我将SVM的所有其他参数保持在默认值。值得注意的是，这意味着惩罚将使用L2准则。

表1详细列出了其余超参数的范围和最终值：标准化方法、特征数量和正则化参数。当然，超参数搜索只在训练集上进行。我使用了五次重复的五倍交叉验证法来估计一组超参数的性能。性能指标是Matthews的相关系数（MCC）。总共尝试了300个组合，MCC为0.478，带有StandardScaler、70个特征、正则化参数C为1.0的配置表现最好。

表1 超参数搜索空间和最终值

虽然我们的主要目标是确定最佳的超参数集，但我们也可以从搜索结果中了解一些关于模型和数据的信息。例如，图3显示了在训练集和测试集上给定固定数量特征的最佳超参数集的性能。虽然数据集的样本比特征多（N>P），但如果我们允许模型使用许多特征，我们观察到典型的过拟合行为。虽然训练集上的性能稳步上升，但测试集上的性能起初缓慢上升，然后趋于平稳，最后随着特征数量的增加而变得更糟。因此，没有特征选择步骤的模型会简单地使用所有可用的特征（图中未显示），其在保持数据上的性能要差很多（MCC 0.414）。

图3 在固定数量的特征下，模型在训练和测试集上的最佳表现

阴影区域代表95%的置信区间，星形表示最佳的整体集合

其他两个超参数，标准化器和正则化参数C，恰好都很不显眼，事实上都处于默认值。

5 模型验证

为了估计模型在生产环境中的真实性能，我们在完整的训练集上使用最佳的超参数训练模型。接下来，我们在一开始就保留的1200个验证样本上评估其性能。在验证集上，SVM管道实现了0.746的准确性和0.480的MCC。这甚至比DNN公布的结果略好。

然而，经过仔细检查，我们发现验证集的组成与最初的问题陈述有冲突。我们的目标是对尚未见过的转录组特征进行分类，但验证集包含许多来自化合物的特征，而这些化合物也是训练集的一部分。因此，我们很可能高估了该模型的真实性能。一个好的经验法则是在交叉验证方案中尽可能地模仿真实的数据采集过程。

在药物发现环境中，这意味着我们有一个具有已知风险特征的化合物池，我们可以用它来训练。然而，该模型随后将被用来估计那些肯定不在训练集中的化合物的DILI风险，因为还没有用这种化合物进行研究或临床试验。因此，交叉验证策略不仅要考虑类的平衡（分层），还要按化合物来拆分数据集。当然，这种逻辑也延伸到了超参数搜索的内部交叉验证。

因此，研究重复了完整的例子，对验证集采用基于化合物的分割，对超参数搜索采用五倍交叉验证的分割方式。不幸的是，在验证集上，准确性下降到0.280，MCC下降到0.113。这表明该模型对未见过的化合物没有很好的概括性，之前的结果是由于该模型记住了单个化合物的特征而取得的。为了检查这种行为是否也发生在原始的DNN模型上，研究对它进行了相同的分割训练。为了得到一个训练好的模型，研究不得不禁用早期停止，因为验证集上的F1分数总是不确定的。因此，最终模型的性能非常低。它的准确率为0.376，MCC为-0.031。总之，这两个模型都不能正确地对新化合物进行分类，因此在真正的药物筛选环境中不会有用。

如果一个项目到了这个地步，我建议退一步，重新考虑一些事情。数据集是否适合解决最初的问题陈述？是否有可能获得更多或更好的数据？其他建模方法是否更有前途？项目的范围是否需要重新定义？也许最好的办法是取消这个项目，重新开始，研究其他的问题。我之所以强调这一点，是因为在组学数据上的机器学习项目有失败的记录，并造成了巨大的声誉和经济损失。我相信，如上所述，通过适当的模型验证，可以大大降低这种失败的风险。

6 最终模型的训练和解释

根据彻底的交叉验证得到的指标，我们可以决定该模型是否满足了继续前进的要求。这可能是一个商业案例的期望，也可能是科学界对学术出版物的期望。为论证起见，我们假设我们对模型的指标感到满意。接下来的步骤是训练最终的模型，了解是什么驱动了模型的分类。幸运的是，第一步是相当容易的。为了生成一个可以在生产中使用的模型，我们需要使用超参数搜索找到的最佳参数最后一次训练管道。这一次，我们使用所有的数据，包括训练、测试和验证集。现在，我们可以使用最终的模型对新的数据集进行预测，如果底层数据源（如测序技术）没有变化，我们可以期待与我们在模型验证期间看到的性能相似。

最后，我们看一下模型的内部运作，弄清楚哪些特征与分类任务有关。做到这一点的一个方法是，从线性SVM中提取每个特征的系数。由于该管道包括一个标准化步骤，所有的特征都在同一个数量级内，因此系数是可以比较的。理论上，我们甚至可以通过系数的符号来确定特征驱动分类的类别。然而，由于特征值可能是负的，解释符号就比较麻烦了，因此我在图4中用绝对系数来显示整体相关性。图中显示了在不按化合物分层的情况下，重复五次的五倍交叉验证中，具有最高中值系数的五个基因的绝对系数分布。

图4 显示五次重复的五倍交叉验证中前五个基因的系数绝对值的博弈图

理想情况下，我们希望看到排名靠前的基因的功能有一个一致的模式。不幸的是，这里的情况并非如此。因此，我编了一个例子。鉴于数据集的性质，基因本体论术语"细胞对药物的反应"的富集本应表明该模型成功地捕获了潜在的生物学意义。然而，我们必须记住，我们正在寻找预测性基因，这些基因不一定对观察到的表型有因果关系。此外，相关特征的系数的权重通常分布在它们之间。显然，这是一个额外的论据，作为特征选择的一部分，聚集或删除相关的特征。

总之，我们看到了数据探索如何指导建模决策，模型定义的哪些方面是重要的，模型验证的陷阱在哪里，以及如何训练和分析最终基于组学的模型。

参考资料

https://link.springer.80599.net/protocol/10.1007/978-1-0716-1787-8_18

----------- End -----------

感兴趣的读者，可以添加小邦微信（zhiyaobang2020）加入读者实名讨论微信群。添加时请主动注明姓名-企业-职位/岗位 或

姓名-学校-职务/研究方向。

- 历史文章推荐 -

【药物设计】